Counting statistics: a Feynman-Kac perspective 
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By building upon a Feynman-Kac formalism, we assess the distribution of the number of hits in 
a given region for a broad class of discrete-time random walks with scattering and absorption. We 
derive the evolution equation for the generating function of the number of hits, and complete our 
analysis by examining the moments of the distribution, and their relation to the walker equilibrium 
density. Some significant applications are discussed in detail: in particular, we revisit the gambler's 
ruin problem and generalize to random walks with absorption the arcsine law for the number of hits 
on the half-line. 



I. INTRODUCTION 

Many physical processes can be described in terms of 
a random walker evolving in the phase space and 
one is often interested in assessing the portion of time ty 
that the system spends in a given region V of the explored 
space, when observed up to time t |5l-[l3t. This is key to 
understanding the dynamics of radiation transport, gas 
flows, research strategies or chemical/biological species 
migration in living bodies, only to mention a few, the 
time spent in V being proportional to the interaction of 
particle with the target medium [l^ - [20j . 

For Brownian motion, a celebrated approach to char- 
acterizing the probability density of the residence time ty 
has been provided by Kac (based on Feynman path in- 
tegrals) in a series of seminal papers, and later extended 
to Markov continuous-time processes [2l| - [2^ . For a re- 
view, see for instance |25i]. Feynman-Kac formalism ba- 
sically allows writing down the evolution equation for 
the moment generating function of ty for arbitrary do- 
mains, initial conditions and displacement kernels. This 
app roach has recently attracted a renovated interest [26l - 
l33j . and has been also extended to non-Markovian pro- 
cesses [13, nil • As a particular case, imposing leak- 
age boundary conditions leads to the formulation of first- 
passage problems 5,, 6,, However, for those physical 
systems that are intrinsically discrete, the natural vari- 
able is the number of hits ny in V when the process is 
observed up to the n-th step, rather than time ty 33, 37|- 
|40| . When ny is large, we can approximate the number 
of events falling in V by ny (x ty (the so-called diffusion 
limit), but this simple proportionality breaks down when 
V is small with respect to the typical step size of the 
walker, and/or the effects of absorption are not negligi- 
ble, so that the diffusion limit is not attained [H, |4l| . 

In this paper, we derive a discrete Feynman-Kac equa- 
tion for the evolution of the probability generating func- 
tion of ny for a broad class of stochastic processes with 
scattering and absorption, and we illustrate this ap- 
proach by explicitly working out calculations for some 
significant examples, such as the gambler's ruin problem, 



or the arcsine law. For the arcsine law, in particular, 
the Feynman-Kac formulae allow generalizing the well- 
known Sparre Andersen results to random walks with 
absorption. Our analysis of the counting statistics is then 
completed by examining the moments of ny, which can 
also be obtained by building upon the Feynman-Kac for- 
malism, and their asymptotic behavior when n is large. 
In particular, we show that the asymptotic moments can 
be expressed as a function of the particle equilibrium dis- 
tribution, which generalizes analogous results previously 
derived in terms of survival probabilities [32] . 

This paper is organized as follows. In Sec. |TT] we in- 
troduce a discrete Feynman-Kac formula for a class of 
random walkers with scattering and absorption. Then, 
in Sec.|nT]we discuss some applications where the gener- 
ating function can be explicity inverted to give the prob- 
ability of the number of hits. In Sec. IIVI we extend our 
analysis to the moments of ny, and in Sec.|V]we examine 
some examples of moment formulae. A short digression 
on the diffusion limit is given in Sec. lVIl Perspectives are 
finally discussed in Sec. IVIII 



II. FEYNMAN-KAC EQUATIONS 

Consider the random walk of a particle starting from 
an isotropic point-source iS(r|ro) = S{r — Tq) located at 
rg. At each collision, the particle can be either scat- 
tered with probability ps, or absorbed with probability 
Pa — 1 — Ps (in which case the trajectory terminates). 
We introduce the quantity T(r' r), namely, the prob- 
ability density of performin g a displacement from r' to r, 
between any two collisions ji^, l^j] • For the sake of sim- 
plicity, we assume that scattering is isotropic, and that 
displacements are equally distributed. 

Suppose that a particle emitted from Fq is observed up 
to entering the n-th collision. Our aim is to characterize 
the distribution P„(nv|ro), where ny is the number of 
collisions in a domain V. We can formally define 



nvin) = ^y(rfc). 



(1) 
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where V{rk) is the marker function of the region V, which 
takes the value 1 when the point rk € V, and vanishes 
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elsewhere. We adopt here the convention that the source 
is not counted, i.e., the sum begins at fc = 1. Clearly, 
ny is a stochastic variable, depending on the realizations 
of the underlying process, and on the initial condition 
Tq. The behavior of its distribution, P„(nv|ro), is most 
easily described in terms of the associated probability 
generating function 



Fn{u\ro) 



■)»(ro) 
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P„(ny|ro)u"^ 



(2) 



which can be interpreted as the discrete Laplace trans- 
form (the transformed variable being u) of the collision 
number distribution. The derivation of an evolution 
equation for ^^(ulro) is made simpler if we initially con- 
sider trajectories starting with a particle entering its first 
collision at ri. Random flights are semi-Markovian (i.e., 
Markovian at collision points) , which allows splitting the 
trajectory into a first jump, from ri to ri -f A (the dis- 
placement A obeying the jump length density T), and 
then a path from ri -|- A to r„, conditioned to the fact 
that the collision at ri was a scattering event. If the col- 
lision is an absorption, the trajectory ends and there will 
be no further events contributing to ny- This leads to 
the following equation for the generating function 



Fn+i{u\ri) = 



,V(ri) 



Ps{F„{u\ri + A)) + Pa] , (3) 



where expectation is taken with respect to the random 
displacement A. The tilde is used to recall that we are 
considering trajectories starting with a single particle en- 
tering the first collision at ri. We make then use of the 
discrete Dynkin's formula 



(/(ri+A)> 



J r*(r' 



ri)/(r')dr'. 



(4) 



where / is any sufficiently well behaved function and 
T*(r' — >■ r) is the adjoint kernel associated to T(r' — > 
r) |4J|. We therefore obtain the discrete Feynman-Kac 
equation 




F„+i(u|ri) 



,V(ri) 



/ T*{r' -^ri)Fn{u\r')dr' +pa 



,(5) 



with the initial condition i^i(M|ri) — u^'^^K Finally, by 
observing that the first collision coordinates ri obey the 
probability density T(ro — ri), it follows 



Fn{u\ro) 



J F„(u|ri)r(ro ^ ri)dri 



FIG. 1: (Color online). The probability P„(Oj3;o) that the 
gambler is not ruined at the n-th bet, given an initial capital 
xo- Bets are modeled by discrete random increments of fixed 
size s = ±1. Blue circles: xo = 5; red triangles: xo = 10; 
green dots: xo = 15. Lines have been added to guide the eye. 
Dashed lines correspond to the asymptotic result P„(0|2;o) — 
s/2/TTxon~^^^ . The interval 2xQ/n is also shown for each xo- 



III. COLLISION NUMBER DISTRIBUTION: 
EXAMPLES OF CALCULATIONS 

Direct calculations based on the discrete Feynman-Kac 
formulae, Eqs. ^ and ([6|), are in some cases amenable 
to exact results concerning P„(ny|ro), at least for simple 
geometries and displacement kernels. In this Section, we 
shall discuss some relevant examples. 



A. The gambler's ruin 

Consider a gambler whose initial amount of money 
is Xo > 0. At each (discrete) time step, the gam- 
bler either wins or loses a fair bet, and his capital in- 
creases or decreases, respectively, by some fixed quan- 
tity s with equal probability. One might be interested 
in knowing the probability that the gambler is not ru- 
ined (i.e., that his capital has not reached zero, yet) 
at the n-th bet, starting from the initial capital xq. 
This well-known problem 



45| can be easily recast in 



(6) 



Knowledge of F„(w|ro) allows explicitly determining 
^n('T-v|ro). Indeed, by construction the probability gen- 
erating function _F„(u|ro) is a polynomial in the variable 
u, the coefficient of each power u'"' being P„(ny = fc|ro). 
In particular, it follows that the probability that the 
particles never touch (or come back to, if the source 
ro e V) the domain V is obtained by evaluating F„(u|ro) 
at u = 0, i.e., P„(0|ro) =F„(0|ro). 



the Feynman-Kac formalism by setting a particle in mo- 
tion on a straight line, starting from xq, with scattering 
probability Ps = I and a discrete displacement kernel 
T{x' x) = S{x - x' - s)/2 + S{x - x' + s)/2. Setting 
s = 1 amounts to expressing the capital xo in multiple 
units of the bet, and entails no loss of generality. The 
counting condition is imposed by assuming a Kronecker 
delta V{x) — 5x.f) in Eq. ([5]): since the walker can not 
cross X = Q without touching it, solving the resulting 
equation for the quantity P'„(0|a;o) gives therefore the re- 
quired probability that the gambler is not ruined at the 
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FIG. 2: (Color online). The arcsine law Pn{nv\xo) with dis- 
crete jump lengths, as a function of nv- The starting point 
is xo = 0, and n = 10. Blue dots: ps = 1 (dots); red stars: 
Ps = 3/4; green circles: Ps = 1/2; black triangles: ps = 1/4. 
Lines have been added to guide the eye. Dashed curves are 
the asymptotic Eq. (|16p for the correponding value of ps- 
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FIG. 3: (Color online). The arcsine law Pn{nv\xo) with dis- 
crete jump lengths, as a function of nv- The starting point 
is Xo = 0, and n = 50. Blue dots: ps — 1; dashed line: 
the asymptotic distribution nv{n — nv)n- Red triangles: 
Ps = 0.95; dashed asymptotic Eq. (fTO)) . 



n-th bet. Now, we integrate Eq. ([SJ and use Eq. ([SJ: 
identifying the coefficient of the zeroth order term in the 
polynomial yields then 

P (n\n\ - /i i ^ ^ A 1£ ^ \ 

' 1 '2'4'8' 16' 32' 64' 128'-/ ^ ' 

for n > 1, when xq = 0. This is easily recognized as 
being the series 



P„(0|0) 



1 



(8) 



[•] denoting the integer part. For xq — 1, 2, 3, ... we ob- 
tain for the first terms in the series 



P„(0|1) = 

Pn(0|2) = 

P„(0|3): 

^n(0|4) 



1 2 3 6 10 20 35 70 
2'4' 8' 16' 32' 64' 128' 256 

2 3 6 10 20 35 70 126 
2'4' 8' 16' 32' 64' 128' 256 
2 4 7 14 25 50 91 182 
2'4' 8' 16' 32' 64' 128' 256 
2 4 8 15 30 56 112 210 
2' 4' 8' 16' 32' 64' 128' 256 



T^'TTTt:,- (9) 



After some rather lengthy algebra, by induction one can 
finally recognize the formula 



Pn(0|a:o) 



r(n+3:o-l)/2] 

E 

k=0 



n 

k — Xq 



2-". (10) 



The quantity Pn{0\xo) is displayed in Fig. [T] as a func- 
tion of n for a few values of xq- The larger the initial 



capital Xo, the longer P„(0|a;o) — 1 before decreasing. 
At large n, taking the limit of Eq. (|10p leads to the 
scaling P„(0|a;o) — y^2/nxon^^^'^ , m agreement with 
the findings in [51]. This means that asymptotically 
the gambler is almost sure not to be ruined, yet, up 
to n ~ 2x'^/tt bets. Note that Eq. ([TOt is the survival 
probability of the gambler: the first passage probabil- 
ity W,i(0|a;o), i.e., the probability that the gambler is 
ruined exactly at the n-th bet, can be obtained from 
Wn{xo) ~ P„_i(0|a;o) — P„(0|a;o). As a particular case, 
for Q < xq < n and n -\- xo even, we recover the result 
in [il], namely 



Wr,{xo) 



(11) 



B. The arcsine law with discrete jumps 

Consider a walker on a straight line, starting from xq- 
We are interested in assessing the distribution Pn{nv\xo) 
of the number of collisions ny that the walker performs at 
the right of the starting point, when observed up to the n- 
th collision. This amounts to choosing V{x) = H{x~xo), 
H being the Heaviside step function, in Eq. ([5]). To fix 
the ideas, without loss of generality we set = 0, and 
we initially assume that = 1, i.e., the walker can 
not be absorbed along the trajectory. This is a well 
known and long studied problern, for both Markovian 
and non-Markovian processes 0,111, m-Sil: for Br own- 
ian motion, the average residence time in V is simply 
(iv)t = whereas ty itself is known to obey the 

so called Levy's arcsine law Pt{ty) = Xj ^ty(t — ty)-K, 
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FIG. 4: (Color online). The arcsine law P„{nv\xo) with ex- 
ponential jump lengths, as a function of ny- The starting 
point is xo ~ 0, and n = 10. Blue dots: Ps ~ 1; red stars: 
Ps = 3/4; green circles: Ps ~ 1/2; black triangles: ps = 1/4. 
Lines have been added to guide the eye. Dashed curves are 
the asymptotic Eq. (|22p for the correponding value of ps- 



FIG. 5: (Color online). The arcsine law P„(nv|a;o) with ex- 
ponential jump lengths, as a function of ny- The starting 
point is Xo = 0, and n = 50. Blue dots: ps = 1; dashed line: 
the asymptotic distribution 1/ ny{n — ny)n. Red triangles: 
Ps = 0.95; dashed line: asymptotic Eq. (|^ . 



whose U shape basically implies that the particle will 
most often spend its time being always either on the pos- 
itive or negative side of the axis [13, IMl . This coun- 
terintuitive result has been shown to asymptotically hold 
also for discrete-time random walks without absorption, 
for which one has P„(ny |0) ~ 1/ y/ny{n — ny)7r when n 
and ny are large (see for instance 

The Feynman-Kac approach allows explicitly deriving 
Pninv\0). Again, assume a displacement kernel with dis- 
crete jumps T{x' — >■ x) = 6{x — x' — s)/2 + S{x — x' + s)/2, 
with s = 1. Then, by integrating Eq. ([S]) and subse- 
quently using Eq. ([6]) we compute the coefficients of 
the polynomial, which can be organized in an infinite 
triangle, whose first terms read 
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Observe that this result is independent of s. To identify 
the elements P„(ny|0), we initially inspect the column 
nv = 1 of the triangle, and recognize the underlying 
series as being given by terms of the kind (|-fe/2])^~'^' 
Then we realize that columns with ny > 2 are related to 
the first column by a shift in the index k. The column 



ny — can be obtained from normalization. Proceding 
by induction, the elements in the triangle can be finally 
recast in the compact formula 



Pniny\0) 



ny — 1 



ny 



(12) 



Note that our result is slightly different from [45|, where 
collisions are counted in pairs. When both n and 
ny are large, we obtain the limit curve P„(ny |0) ~ 
^/\/ny{n — ny)TT. 

When the scattering probability can vary in 
< < 1, the triangle can be generated as above, and 
the first few terms read 
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Now, the identification of the polynomial coefficients 
P„(ny|0) becomes more involved, because each coeffi- 
cient is itself a polynomial with respect to Ps ■ The strat- 
egy in the identification is the same as above. By induc- 
tion, the column for ny — 1 can be identified as being 



Pnim 



p x4+2W2 + 2y\ 2i^l (l,i+y,3 + y,p2) 



(I) 



2(2 + 1/) 



(13) 



for n > 2, with y = \{n - 3)/2], and Pi(l|0) = 1/2. 
Once P„(1|0) is known, by inspection one realizes that 
the other columns Pn{nv\^) are related by 



Ps = 3/4 



Ps = 3/4 



Pn{nv\(^)^{^) 



nv 



Pn- 



riv + l 



(1|0) (14) 



for nv > 2. The probability P„(0|0) is finally obtained 
from normalization, and reads 



P„(0|0) 



Ps-1 + v/l-Ps 



2+2z /2z 

z 



z,l,2 + z,p2) 



p,(l + z) 



(15) 



with z — \nl2\. These results generalize Eq. ([T^. and 
are illustrated in Fig. [2j where we compare P„(n\/|0) as 
a function of ny for n = 10 and different values of p^. 
When Ps = \ the distribution approaches a U shape, as 
expected. As soon as < 1, the shape changes con- 
siderably, and in particular P„(ny|0) becomes strongly 
peaked at ny = as the effects of absorption overcome 
scattering. The presence of a second peak at ny = n 
is visible when ~ 1 and progressively disappear as 
decreases: whenps is small, P„(nv'|0) has an exponential 
tail. When n is large, P„(nv|0) approaches the asymp- 
totic curve 



i^oo(nv|0)=(|)' 



ny 



1 - p. + Vi-p' 



(16) 

for ny > 1, and Poo(0|0) = (p. - 1 + v/T^)/2p.. 
Remark that when ps = 1 this means that the U shape 
for large n collapses on the two extremes at ny = and 
ny = n. Equation ()16|) is an excellent approximation of 
Pn{ny\Q) when the scattering probability is not too close 
to Ps — I: as expected, the discrepancy between the 
exact and aymptotic probability is most evident when 
ny ~ n, as shown in Fig. [2] Fig. [3] displays P„(nv|0) as 
a function of ny for n = 50, in order to emphasize the 
effects of p^: when ps = 1 the probability P„(ny|0) is 
almost superposed to the asymptotic curve P„(ny|0) ~ 
1/ y/ny{n — ny)Tr, whereas a deviation in the scattering 
probability as small as Ps = 0.95 is sufficient to radically 
change the shape of the collision number distribution. 
Finally, observe that when ny is also large, which implies 
Pa <^l, Eq. (|T6)) behaves as 



Poo(ny|0) 



nny 



(l-Ps)nv 



(17) 



C. The arcsine law with continuous jumps 

When the displacement kernel T{x' — >■ x) is continuous 
and symmetric, = 1, and xq — 0, the distribution of 
the number of collisions ny falling in x > xq is universal. 
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FIG. 6: (Color online). One-dimensional exponential flights 
in V — [—1, 1]. Upper half. First (left) and second factorial 
moment (right) of nv as a function of n and xq, with ps =3/4. 
The asymptotic values for a;o = are shown as a dahsed lines 
for reference. Lower half. First (left) and second factorial 
moment (right) of ny in the same domain, when ps = 1. 



in that it does not depend on the specific functional 
form of T{x' — >■ x) (see [5l'| and references therein). 
This strong and surprising result stems from a Sparre 
Andersen theorem |52l | , whose proof is highly non trivial 
(and does not apply to discrete jumps) [51|. This leaves 
the choice on the form of kernel T{x' — )■ a;), as far as it 
satisfies the hypotheses of the theorem. For the sake of 
simplicity, we have assumed an exponential distribution 
of jump lengths, i.e., T(x' x) ~ sexp(— sja: — x'|)/2, 
with s — \. Starting from Eq. ([5]) and Eq. ^ we can 
generate the infinite triangle 
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It is easy to verify that the triangle indeed does not 
depend on s, and that other functional forms of T(x' 
x) would lead to the same coefficients for the polynomials 
P„(n\/|G). This holds true also for Levy flights, where 
T{x' x) is a Levy stable law and jump lengths are 
unbounded [5l|. We start from the column ny = 1, 
observe the relation with the subsequent columns ny > 
2, and finally derive the case ny = from normalization. 
Proceding therefore by induction we recognize that the 
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elements of the triangle obey 

P„(n.|0)=P"~'"^)P"^)2--. (18) 
V n-nv J \nv J 

We recover here the celebrated results of the collision 
number distribution for discrete-time walks with sym- 
metric continuous jumps, in absence of absorption [i^, 
[Slj . When both n and ny are large, it is possi- 
ble to show that Eq. (fTSi) converges to the U shape 
1/ y^riy(n^^rT^7r. 

When the scattering probability is allowed to vary in 
< Ps < 1, it turns out that the polynomial coefficients 



P„(ny|0) are the same for several different continuous 
symmetric kernels T{x' — >■ x) (Levy flights included), 
and we are therefore led to conjecture that the uni- 
versality result for the case ps — 1 carries over to to 
random walks with absorption. This allows generalizing 
the Sparre Andersen theorem for the collision number 
distribution on the half-line to a broader class of Marko- 
vian discrete-time processes. The first few terms in the 
triangle (which for practical purposes we have generated 
by resorting to T{x' x) = sexp{—s\x — x'\)/2, with 
s — 1) read 
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As above, identification of the terms P„(ny |0) becomes 
more involved, because each coefficient is itself a polyno- 
mial with respect to Ps- By induction, the column for 
nv = 1 can be identified as being 



Pnim 



Ps\"/2n-2\ 2Fi{-^+n,l,l + n,ps 



n — 1 



(19) 



for 71 > 1. Once P„(1|0) is known, the subsequent 
columns Pn{nv\0) are observed to obey 

PninvlO) = (^"l^ ^)p„-„,,+i(l|0) (20) 

for nv > 2. The probability P„(0|0) is finally obtained 
from normalization, and reads 



P„(0|0) 



Ps 



1 + 



Ps 



2n\ 2F1 {} + n,l,2 + n,Ps 
2(1 + n) 



(21) 



These results are illustrated in Fig. 21 where we compare 
Pri("-y|0) as a function of nv for n = 10 and different 
values of Ps- The findings for continuous jumps closely 
resemble those for discrete displacements. When ps — 1 
the distribution approaches a U shape, as expected. As 
soon as Ps < 1, the shape again changes abruptly, and in 
particular P„(ny |0) becomes strongly peaked at nv — 
when absorption dominates scattering. When ps is small. 



P„(ny|0) decreases exponentially at large ny. When n 
is large, P„(ny|0) approaches the asymptotic curve 



Poo(«l^|0) = 



nv 



2nv - 1\ \/l - p. 



(22) 



forny > 1, andPoo(0|0) = (ps-l + \/l -Ps)/Ps- Again, 
Eq. is an excellent approximation of P„(ny|0) when 
the scattering probability is not too close to Ps ~ 1, as 
shown in Fig. SI Fig. [S] displays P„(ny |0) as a func- 
tion of nv for n = 50, in order to emphasize the ef- 
fects of Ps'- when ps = 1 the probability P„(ny |0) is 
almost superposed to the asymptotic curve P„(ny |0) ~ 
1/ y^nv{n — nv)iT, whereas a deviation in the scattering 
probability as small as Ps — 0.95 is sufficient to radically 
change the shape of the collision number distribution. Fi- 
nally, observe that when ny is also large, which implies 
Pa < 1, Eq. dm yields the same scaling as Eq. ([T71). All 
analytical calculations discussed here have been verified 
by comparison with Monte Carlo simulations with IO6 
particles. 



IV. MOMENTS FORMULAE 

A complementary tool for characterizing the distribu- 
tion P„(ny |ro) is provided by the analysis of its moments. 
To this aim, it is convenient to introduce the function 
G„('u|ri) — P„(l/M|ri). By construction, G„(u|ri) is 
the (rising) factorial moment generating function for tra- 
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jectories entering their first collision at ri, which implies 



(nr)«(ri) = (-1)' 



rG„(u|ri)U=i, (23) 



x{x + 1) ...{x + k — 1) being the rising factorial [52 



The tilde is used to recall that the moments refer to tra- 
jectories starting from the first collision at ri. Combining 
Eqs. (O and ([25)1 yields the recursion property 



r,(™)\ 



V 



)n+l 



(ri),(24) 



for 



> 1, with the conditions 



(ri) — 1, and 



(nl7*'')i(ri) = m\V(ri). Finally, the factorial moments 



(ny"^)„(ro) for particles emitted at ro are obtained from 



/n 



(ro) 



-(™)\ 

/n 



(ri)T(ro ^ ri)dri. (25) 



When trajectories are observed up to n — > +oo, we can 
set (riy"^) = lim„_j.+oo(»^v™^)n, and from Eq. (p4)) we find 



(nlr^)(ri) -P.y T*(r' -> ri)(ni;'^>(r')dr' = 

= mF(ri)(n(r"^>)(ri), (26) 

for m > 1, provided that {h^^){ri) is finite. It turns 

out that the asymptotic moments (n(7''')(ri) are related 
to the equilibrium distribution of the particles [32l . Issl . 
Isgj . To see this, we introduce the incident propagator 
'J'„(r|ro), i.e., the probability density of finding a particle 
emitted at Tq entering the n-th collision (n > 1) at r. We 
have 



*„+i(r|ro) =PsJ T(r' ^ r)^'„(r'|ro)dr'. 



(27) 



with ^'i(r|ro) = T(ro r). We introduce then the inci- 
dent collision density 



TV 



*(r|ro) = Jim V *„(r|ro), 

A^— 5-00 ^ — ' 



(28) 



n=l 



which can be interpreted as the particle equilibrium dis- 
tribution [s^, m, 112, E3| • Now, by making use of the for- 
mal Neumann series [43| , from Eq. (P7|) it follows that the 
collision density satisfies the stationary integral transport 
equation 

^'(r|ro) =psj T{r' ^ r)*(r'|ro)rfr' + r(ro ^ r). (29) 



The transport equation p9| can be understood as fol- 
lows: at equilibrium, the stationary particle density en- 
tering a collision at r for a source emitting at Tq is given 
by the sum of all contributions entering a collision at r'. 



being scattered and then transported to r, plus the con- 
tribution of the particles emitted from the source and 
never collided up to entering r. Now, by resorting to the 
relation between T and its adjoint T* in Eq. (|A4[) . and 
observing that at equilibrium (for isotropic source and 
scattering) 

J r(ro ^ r')*(r|r')dr' = J *(r'|ro)r(r' ^ r)dr', 

(30) 

then Eq. (|26p can be inverted (see Appendix and 
gives 



Jv 

+mV{r,){n'(^-'^){r,), (31) 



for TO > 1, where ^E* is the solution of Eq. ([29)) . Thus, by 
using Eqs. (US]) and dSHl), from Eq. it follows 



(4™))(ro) = m/ vl,(r'|r„)(n;;'-^^)(r')dr'. (32) 
Jv 

Hence, by induction we finally get the desired relation 
between the factorial moments and the equilibrium dis- 
tribution, namely. 



^'lr.„^/^^("'-l)\('v.'^Wr.' 



)(ro) =Y,Lm,kPs ^Cfc(ro), 



fc=i 



where 



Ck{ro) = k\ 



J|dr,*(r,|r,_i) 

i=l 



(33) 



(34) 



are fc-fold convolution (Kac) integrals over the equilib- 
rium distribution \I'(r|ro), and the coefficients Lm.fc = 
(")(m - l)!/(fc - 1)! are the Lah numbers jH. 'Ob- 
serve that a result analogous to Eq. ([55)1 has been derived 
in [32] building upon survival probabilities. 



V. FIRST AND SECOND MOMENT OF 
COLLISION NUMBER: EXAMPLES OF 
CALCULATIONS 

To illustrate the approach proposed in the previ- 
ous Section, we compute here the average collision 

number (ny'')„(a;o) and the second factorial moment 

(riy )„(a;o) in a bounded domain V in one dimension. 
We choose an exponential displacement kernel T{x' — > 
x) — sexp(— s|a; — x'|)/2, with s = 1, which is often 
adopted as a simplified random walk model (the so-called 
exponential fiights) to describe gas dynamics, radiation 
propagation or biological species migration |54| - [59| . In 
this respect, the average (ny-')„(xo) is a measure of the 
passage of the particles through the region V (the de- 
posited energy, for instance) , whereas the second moment 
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.(2) 



)n(xp) is proportional to the incertitude on the av- namely, 



erage [32|, |33[. 

The calculations stemming from Eqs. (1^^ and (pSj) are 
rather cumbersome, so that it is preferable to visually 
represent our results instead of writing down the explicit 
formulas. In Fig. [6] we display the behavior of the mo- 



ments {ny^)nixo) and {ny')nixa) for the volume V being 
the interval [—R, R], with R — 1. When ps = 1, the mo- 
ments diverge as n increases, since exponential flights in 
one dimension are recurrent random walks, and revisit 
their starting point infinitely many times. When < 1, 
they converge instead to an asymptotic value, which can 
be computed based on Eq. (l33t by observing that the 
collision density for this example is 



,(2) 



«'(x|a;o) = 



-^/l-ps\x-Xo\ 



2VT 



Ps 



(35) 



as discussed in [32]. Moreover, Fig. |6] shows that the 
moments decrease as the distance of the source Xo from 
the region V increases, as expected. As we have cho- 



sen here a symmetric interval, we have 

(m 



nixo) 



(riy )n{—xo), so that we can plot the moments only for 
positive values of xq. All results presented in this Sec- 
tion have been verified by comparison with Monte Carlo 
simulations with 10^ particles. Other kinds of boundary 
conditions (leakage, for instance, which implies that par- 
ticles are lost upon crossing the frontier of V [H, HI]) 
have also been successfully tested, but will not be pre- 
sented here. 



VI. DIFFUSION LIMIT 

To conclude our analysis, in this Section we comment 
on the scaling limit of the discrete Feynman-Kac equa- 
tion, which is achieved when ny is large, and at the same 
time the typical jump length e is vanishing small. We set 
tv = nydt and t = ndt, where dt is some small time scale, 
related to e by the usual diffusion scaling ~ 2Ddt, 
the constant D playing the role of a diffusion coefficient. 
When T is not symmetric, so that displacements have 
mean /x, we further require /x = vdt, where the constant 
u is a velocity. By properly taking the limit of large ny 
and vanishing dt, ty converge to the residence time in 
V . The quantity ny can only be large if the absorption 
probability Pa is small, and it is natural to set Pa = Xadt, 
the quantity Xa being an absorption rate per unit of dt. 
Observe that when both e and /i are small for any dis- 
placement kernel we have the Taylor expansion 

J T(v' ^ r)/(r')dr' ^ /(r) - ^d.f{v) + \e^dlf{v), 

(36) 

where the first-order derivative vanishes if the kernel is 
symmetric. A similar expansion holds for the kernel T*, 



r*(r' ^ ro)/(r')dr' ^ f{vo)+^,d.J{v^)+-e'dlj{vo). 

(37) 

It is expedient to introduce the quantity (5t(u]ro) = 
F4(e~"]ro), which is the moment generating function of 
ty = nydt, i.e.. 



(C)*(ro) = (-l)^ 



Qr, 



Mu\ro)\u=o, (38) 



when trajectories are observed up to t = ndt. Under the 
previous hypotheses, combining Eqs. ([5]) and ^ yields 

Qt+dt{u\ro) - Qt{u\ro) ~ 
~ Cl^Qtiu\ro)dt ~ uV{r„)Qt{u\ro)dt + Kdt, (39) 

where we have neglected all terms vanishing faster than 
dt, and = Dd^^ + vdr^ - Xa- Taking the limit dt —i' 
we recognize then the Feynman-Kac equation for a Brow- 
nian motion with diffusion coefficient D, drift v and ab- 
sorption rate Xa, namely 



dQt{u\ro) 
dt 



C*^^Qt{u\ro) - uU(ro)Qt(u]ro) + Xa. (40) 



In other words, in the diffusion limit the statistical prop- 
erties of the hit number in V behave as those of the res- 
idence time of a Brownian motion, as is quite naturally 
expected on physical grounds [13, [H, [slj . Finally, from 
Eq. (j38p stems the recursion property for the moments 



C:,{f{^) tiro) + mViro){t"^-') tiro), (41) 



^(^y)t(ro) 
dt 



in agreement with the results in [29l l30j for Brownian 
motion. 



VII. CONCLUSIONS 

In this paper we have examined the behavior of the 
distribution P„(ny]ro) of the number of hits ny in a 
region V for a broad class of stochastic processes with 
scattering and absorption. Key to our analysis has been 
a discrete version of the Feynman-Kac formalism. We 
have shown that this approach is amenable to explicit 
formulas for P„(nv]ro), at least for simple geometries and 
displacement kernels. The moments of the distribution 
have also been detailed, and their asymptotic behavior 
for large n has been related to the walker equilibrium 
density. Finally, the diffusion limit and the convergence 
to the Feynman-Kac formulae for Brownian motion have 
been discussed. 

We conclude by observing that a generalization of the 
present work to more realistic transport kernels, includ- 
ing anisotropic source and scattering, would be possi- 
ble, for instance by resorting to the formalism proposed 
in [SQl- Moreover, while in this paper we have focused 
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on counting statistics, and therefore chosen ^(r) to be 
the marker function of a given domain in phase space, 
the Feynman-Kac formahsm can be adapted with mi- 
nor changes to describing the statistics of other kinds 
of functionals, such as for instance hitting probabih- 
ties [i,|3i,|6li. 



Appendix A: The stationary moment equation 

We want to solve an integral equation 

/(n) -psj T*{r' ^ ri)/(r')dr' = g{r,) (Al) 

for the function /(ri), where g{ri) is known. We propose 
a solution in the form 

/(n) - y 5(r')*(r'|ri)dr' + ^(ri) (A2) 

and ask which is the equation satisfied by the integral 
kernel *(r'|ri). By injecting Eq. (|X2|) into Eq. (KJ) one 
obtains 



5(r')*(r'|ri)dr' = 

^{r"\r')T*{r' ri)g{r")dr' dr" + 

+ J T*{r' ^ ri)g{r')dr'. (A3) 

Recall that the adjoint and direct displacement kernels 
are related to each other by the scalar products 

Jg{r) j T{v' ^v)f(v')dv'dv = 

= j /(r) j T*{r' ^r)g{r')dr'dr (A4) 



for any test functions / and g . From the definition of 
the scalar product in Eq. (|A4p it follows that the second 
term at the right hand side of Eq. (jA3|) is given by 



J T*{r' -> ri)g{r')dr' ^ J T{ri ^ r')g{r')dr'. (A5) 



From Eqs. (|A4[) and ([501 the first term at the right hand 
side of Eq. (|A3[) becomes 



*(r"|r')r*(r' ^ ri)g{r")dr' dr" = 

^'(r'lri)T(r' ^ r")ff(r")dr'dr". (A6) 



Therefore, 4'(r'|ri) obeys 

J g(r')*(r'|ri)dr' = 
= Ps f f *(r'|ri)T(r'^r")g(r")rfr'dr" + 



+ J T(ri ^ r')5(r')dr', (A7) 

which for the arbitrariness of 5(r') finally implies 
Eq. ([29l) . i.e., the required kernel \l/ satisfies the integral 
transport equation. 
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